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ABSTRACT':  Zero-crossing  or  feature-point  based  stereo  algorithms  can,  by  definition,  determine 
explicit  depth  information  only  at  particular  points  in  the  image.  To  compute  a  complete  surface 
description,  this  sparse  depth  map  must  be  interpolated.  A  computational  theory  of  this  interpolation 
or  reconstruction  process,  based  on  a  surface  consistency  constraint,  has  previously  been  proposed. 
In  order  to  provide  stronger  boundary  conditions  for  the  interpolation  process,  other  visual  cues 
to  surface  shape  are  examined  in  this  paper.  In  particular,  it  is  shown  that,  in  principle,  shading 
information  from  the  two  views  can  be  used  to  determine  the  orientation  of  the  surface  normal 
along  the  feature-point  contours,  as  well  as  the  parameters  of  the  reflective  properties  of  the 
surface  material.  'Hie  numerical  stability  of  the  resulting  equations  is  also  examined. 
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part  by  the  Advanced  Research  Projects  Agency  of  the  Department  of  Defense  under  Office  of 
Naval  Research  contract  N00014-80-C-0505  and  in  part  by  National  Science  Foundation  Grant 
79-231 10MCS. 


I.  Introduction 


1.1.  The  Global  Problem 

One  of  the  goals  of  a  vision  system  is  to  compute  the  three-dimensional  shape  of  the 
visible  surfaces  in  a  scene.  The  human  visual  system  uses  many  cues  to  compute  surface  shape, 
with  different  modules  of  the  system  using  varying  sources  of  information  in  the  images  to 
infer  information  about  surface  shape.  Hxamplcs  include  motion  analysis,  stercopsis.  shading,  and 
texture.  How  do  these  different  visual  modalities  contribute  to  the  computation  of  surface  shape? 

One  computational  approach  to  understanding  the  human  \  isual  system,  pioneered  by  Marr 
and  Poggio  (sec,  for  example.  [Marr,  1976,  1982:  Marr  and  Poggio,  1977]).  views  the  computation 
of  surface  shape  from  images,  in  part,  as  a  collection  of  transformations  between  two  main 
representations,  Ihe  first  representation  is  the  primal  sketch ,  which  makes  explicit  loci  of  changes 
in  image  irradiance  at  particular  scales  of  resolution:  the  second  is  the  2i -I)  sketch,  which  makes 
explicit  information  about  surface  shape  and  reflective  properties  of  the  surface  material.  The 
modules  that  compute  information  feeding  the  2  £-1)  sketch  from  the  primal  sketch  have  generally 
been  considered  to  a  first  approximation  to  be  independent  of  one  another.  It  is  clear,  however, 
that  within  the  2.J-I)  sketch  the  different  sources  of  information  should  interact,  both  to  maintain 
consistency  among  the  data  provided  by  different  modules,  and  to  provide  feedback  to  the  modules 
in  order  to  enhance  the  acquired  data  (for  example,  texture  contours  can  facilitate  stercopsis 
by  driving  vcrgcncc  eye  movements  [Kidd,  et  al.,  1979]).  In  this  paper,  we  arc  interested  in 
examining  interactions  at  the  level  of  the  2$-D  sketch  between  modules  of  the  early  visual  system. 
In  particular,  wc  will  investigate  some  of  the  ways  in  which  shading  information  can  augment 
stereo  data. 

1.2.  The  Motivating  Problem 

The  goal  of  the  2J-D  sketch  [see,  for  example,  Marr.  1978.  1982]  is  to  compute  surface 
parameters,  in  particular,  the  distance  to  and  orientation  of  small  patches  of  the  visible  surfaces, 
the  discontinuities  in  those  surfaces  (for  example,  the  edges  of  objects),  and  possibly  the  properties 
of  the  surface  material  (for  example  the  amount  of  specularity,  the  colour  and  the  albedo  of 
the  surface  material).  Representations  similar  to  the  2^-D  sketch  have  also  been  suggested  by 
Horn  [1979]  and  Barrow  and  Tcnenbaum  [1979],  As  mentioned,  the  input  to  this  representation  is 
provided  by  a  number  of  roughly  independent  modules,  two  important  ones  of  which  are  stercopsis 
and  motion  perception.  The  input  provided  by  these  particular  modules,  which  is  characteristic 
of  many  early  visual  modules,  consists  of  explicit  information  about  surface  shape  or  disposition 
only  at  particular  points  in  the  image. 

This  follows  from  the  form  of  the  primal  sketch.  Both  psychophysical  investigations  [for 
example,  Attneave,  1954;  Barlow,  1961]  and  computational  studies  [for  example,  Roberts,  1965; 
Hcrskovitz  and  Binford,  1970;  Horn,  1973:  Marr,  1976;  Marr  and  Hildreth,  1980,  and  reviews 
by  Roscnfcld  and  Kak,  1976;  and  Pratt,  1978]  of  early  visual  processing  suggest  that  most  of 
the  information  in  an  image  is  carried  by  the  intensity  changes,  and  hence,  that  the  first  stage 
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of  analysis  of  images  is  the  detection  of  such  changes.  These  changes  in  intensity  arc  used  as 
input  by  the  modules  that  feed  the  2^-1)  sketch  and  hence  explicit  information  is  obtained  only 
at  those  locations.  For  example,  in  stereo  computations  [Marr  and  Foggio,  1979;  Crimson.  1981; 
Baker  and  Binford.  1981;  May  hew  and  Frisby,  1981]  the  zero-crossing  contours  of  the  primal 
sketches  [Marr  and  Hildreth.  1980]  for  the  left  and  right  eye  are  placed  in  correspondence,  and 
the  difference  in  projection  of  the  corresponding  contours  in  the  two  eyes  is  used  to  determine  the 
depth  of  the  surface  along  that  contour.  To  create  a  complete  surface  representation  it  is  necessary 
to  interpolate  between  the  known  values  of  the  raw  2^-D  sketch.  An  initial  theory  of  this  process 
that  implicitly  takes  into  account  some  of  the  shading  information  available  in  the  primal  sketch 
has  been  proposed  [Grimson.  1982a,  1982b;  see  also  refinements  by  Tcr/opoulos,  1982], 

While  this  interpolation  or  reconstruction  theory  constructs  die  “best"  surface  to  fit  through 
the  known  depth  points  provided  by  stereopsis.  it  is  clear  that,  in  prin-'1- supplying  additional 
constraints  to  the  reconstruction  process  would  result  in  an  improved  accuracy  in  the  computed 
surface  shape. 

An  extreme  example  of  this  can  be  seen  in  Figure  1,  which  shows  one  possible  set  of  depth 
constraints  obtained  by  stcrcoscopicaily  viewing  a  sinusoidal  depth  grating  of  a  uniform  material 
from  a  particular  viewpoint.  Ihc  most  conservative  reconstruction,  without  explicitly  taking  the 
shading  information  into  account,  would  be  the  plane  shown  in  Figure  lb.  On  the  other  hand, 
if  surface  orientation  information  were  also  available  along  the  depth  contours,  the  more  correct 
reconstruction  of  Figure  lc  would  be  obtained.  Ihus.  the  motivation  is  to  find  algorithms  for 
computing  surface  orientation  at  the  known  depth  points  prov  ided  by  stereopsis. 

1.3.  The  Specific  Problem 

Given  that  we  want  to  augment  the  boundary  conditions  supplied  to  the  reconstruction 
process,  the  specific  question  to  be  investigated  here  is  whether  shading  information  along  the 
zero-crossing  contours  of  the  primal  sketch  can  be  used  to  provide  information  about  the  surface 
orientation  along  those  contours. 

The  traditional  approach  to  the  shape-from-shading  problem,  pioneered  by  Horn  [1970,  1975], 
has  been  monocular.  The  goal  has  been  to  extract  explicit  surface  orientation  values  at  each  point  in 
the  image,  from  a  single  view  of  the  scene.  As  stated,  the  problem  is  considerably  underconstrained 
and  thus  additional  information  is  required  in  order  to  obtain  the  surface  orientation  information 
(sec  for  example  [Ikcuchi  and  Horn,  1981;  Bars,  1981;  Brooks.  1982]).  For  example,  it  is  usually 
assumed  that  the  direction  and  strength  of  the  illuminant  and  the  reflective  properties  of  the 
surface  material  (which  are  assumed  constant  over  large  sections  of  the  image)  arc  known.  One 
method  for  obtaining  the  additional  information  necessary  to  solve  for  the  surface  orientation, 
under  these  assumptions,  is  to  use  the  technique  of  photometric  stereo  [Woodham,  1978,  1980, 
1981;  Ikcuchi  and  Horn,  1979;  Silver,  1980].  In  this  ease,  multiple  images  obtained  from  a  single 
viewpoint,  but  with  different  illuminant  directions,  are  used  to  determine  the  surface  orientation 
at  points  in  the  image. 

The  intention  in  this  paper  is  to  solve  (at  least  partially)  the  shading  problem  with  fewer 
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Figure  1.  Part  (a)  shows  an  extreme  example  of  a  scl  of  depth  constraints  that  could  be  obtained  by 
stereoscopical  I  y  viewing  a  sinusoidal  depth  grating.  If  shading  information  is  not  explicitly  accounted  for, 
the  best  surface  reconstruction  is  the  plane  shown  in  part  (b).  On  the  other  hand,  if  the  shading  information 
is  used  to  obtain  surface  orientation  information  along  the  depth  contours,  the  more  accurate  reconstruction 
of  part  (c)  is  obtained. 


requirements  or  assumptions  about  the  parameters  of  the  imaging  process.  To  Jo  so  v*ill  require 
either  additional  sources  of  information  about  shading  or  weaker  expectations  on  the  delivered 
data,  or  both.  Indeed,  since  the  motivating  problem  is  to  obtain  additional  boundary  conditions  on 
the  reconstruction  process,  we  will  only  expect  surface  orientation  information  at  the  zero-crossings, 
rather  than  everywhere  in  the  image.  We  will  further  only  expect  a  coarse  estimate  of  the  surface 
normal,  since  the  values  will  serve  as  boundary  conditions  for  a  surface  approximation  process, 
rather  than  constituting  an  explicit  surface  reconstruction  in  itself.  Moreover,  we  will  introduce 
an  additional  source  of  information  to  the  process.  In  particular,  we  will  investigate  the  added 
information  provided  by  the  differing  shading  information  observed  at  nearby  viewpoints.  Hence, 
we  are  concerned  with  binocular  shading,  especially  as  it  applies  to  the  surface  reconstruction 
problem.  The  main  questions  to  be  addressed  here  are:  (1)  Using  binocular  shading,  under  what 
conditions  can  we  solve  for  the  surface  orientation  at  the  zero-crossing  contours  that  have  been 
matched  by  stcreopsis,  and  (2)  How  numerically  stable  is  the  resulting  process?  We  will  see  that 
the  method  used  to  solve  for  the  surface  orientation  is  essentially  an  alternate  fonn  of  photometric 
stereo,  using  images  obtained  from  different  viewpoints,  rather  than  with  different  illuminants. 

2.  The  Basic  Equations 

The  brightness  recorded  at  a  point  in  an  image  can  generally  be  considered  as  a  product 
of  three  factors,  the  total  amount  of  incident  light  striking  the  surface,  the  percentage  of  such 
incident  light  which  is  reflected  (as  opposed  to  absorbed  or  transmitted),  usually  denoted  by  the 
albedo,  p.  of  the  surface,  and  the  reflectance  of  the  surface,  the  distribution  of  the  reflected  light  as 
a  function  of  direction.  Here,  we  shall  assume  that  the  intensity  of  the  incident  light  is  normalized, 
by  incorporating  it  into  the  definition  of  p.  As  a  consequence,  the  parameter  p  is  not  restricted  to 
the  range  0  —  1,  but  rather  can  take  on  any  positive  value. 

Using  a  right-hand  coordinate  system  with  origin  at  the  lefl  eye,  x  axis  connecting  the  two 
eyes  and  the  negative  z  axis  pointing  straight  ahead,  a  surface  may  be  represented  by  a  Monge 
patch 

{*,y, -/(*.»)} 

and  the  surface  normal  at  a  point,  by 

{P.9.1} 

where 

0=*/  o=V 

dx  dy 

The  situation  to  be  investigated  here  is  one  in  which  the  illumination  geometry  consists  of  an 
arbitrary  surface  illuminated  by  a  point  source,  sufficiently  distant  from  the  surface  relative  to  the 
distance  to  the  viewer.  Under  most  circumstances  (for  example,  the  surface  reflectance  is  isotropic), 
the  reflectance  properties  of  the  surface  can  be  combined  with  the  illumination  geometry  and  the 
viewer  position  into  the  convenient  representation  of  a  reflectance  map  [Horn,  1977). 

For  an  arbitrary  surface,  the  surface  reflectance  map  will  generally  be  a  combination  of  two 
effects.  The  first  is  a  diffuse  or  matte  component,  which  results  from  the  scattering  of  light  that 
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penetrates  some  distance  into  the  microstructure  of  the  surface.  The  second  is  a  specular  or  luster 
component,  which  results  from  the  mirror-like  reflection  of  light  at  a  smooth  interface  between 
two  materials  of  different  refractive  index.  These  tw'o  components  can  be  combined  into  a  single 
reflectance  map  in  the  following  manner  (see,  for  example.  (Horn,  1981;  Blinn,  1977;  Phong, 
1975]). 

I  .ct  n  denote  a  unit  vector  in  the  direction  of  the  surface  normal,  defined  at  a  point  on  the 
surface,  let  s  denote  a  unit  vector  in  the  direction  of  the  point  light  source,  and  let  v  denote  a 
unit  view  vector.  (In  general,  we  will  consider  objects  to  be  distant  relative  to  the  viewer,  so  that 
v  is  essentially  constant,  and  will  consider  objects  to  also  be  distant  relative  to  die  light  source,  so 
that  s  is  essentially  constant).  Ihcn  die  recorded  brightness  is  given  by 

c  =  pm 

and  the  reflectance  map  is  given  by 

rt(n)  =  {1  —  a)(ns)  +  a{n  •  h}*  (1) 

where  h  is  a  unit  vector  in  die  direction  of  the  off  specular  angle. 


Here,  p,  a  and  k  arc  parameters  of  the  particular  surface,  with  p  and  k  non-negative  and  a  ranging 
between  0  and  1. 

Note  that  n  is  a  function  of  die  surface  /  while  s  is  a  function  of  the  illumination  geometry. 
Tor  our  particular  choice  of  coordinate  system,  the  view  vector  for  the  left  eye  is  particularly 
simple,  V;.  =  {0,0,1}.  Now  suppose  that  we  take  a  second  view  of  the  surface,  with  vergence 
angle  0  (see  Figure  2).  Here,  the  view  vector  is  given  by  V/*  =  {— sin/3,0,  cos/3}.  Note  that  the 
vergence  angle  0  and  the  illuminant  direction  are  assumed  to  be  known,  so  that  the  vectors  vR 
and  s  arc  constants. 

The  first  term  in  expression  (1)  corresponds  to  a  matte  or  I^mbcrtian  component  of  the 
surface.  Since  this  component  depends  only  on  the  angle  between  the  light  source  and  the  surface 
normal  (and  hence  is  independent  of  v),  it  is  identical  for  both  views.  Thus,  no  additional 
information  will  be  obtained  from  two  views  in  the  ease  of  a  perfectly  matte  surface.  The  second 
term  corresponds  to  a  specular  component  of  the  surface,  and  this  term,  in  general,  does  change 
oetween  the  two  views,  as  can  be  seen  in  expression  (1). 

Our  intent  is  to  use  the  two  expressions  for  recorded  brightnesses  in  the  two  eyes  to  determine 
the  orientation  of  the  viewed  surface  at  some  point.  In  order  to  compare  irradiance  values  from  the 
two  views,  we  must  be  certain  that  the  points  in  the  images  at  which  the  irradiance  measurements 
arc  taken  arc  in  correspondence  (that  is,  we  know  which  points  in  the  two  images  correspond  to 
the  same  point  on  the  surface).  In  the  ideal  ease  of  isolated  intensity  edges,  matched  zero-crossings 
from  a  stereo  algorithm  (for  example,  Grimson,  1981]  are  in  correspondence,  so  by  restricting  our 
attention  to  such  points,  we  arc  in  principle  guaranteed  irradiance  measurements  that  correspond 
to  the  same  point  on  a  physical  surface.  We  should  note,  however,  that  in  practice  a  certain  amount 
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Figure  1  The  geomeiry  of  iwo  views  of  ihe  suiface.  with  the  vergence  angle  0  indicated. 


of  error  will  be  associated  with  the  matched  zero-crossings,  due  to  uncertainly  in  exactly  localizing 
the  underlying  intensity  changes.  This  error  will  reduce  the  effectiveness  of  the  computation. 

There  are  two  stages  to  the  computation  of  surface  orientation,  ITie  first  is  to  compute  the 
surface  reflectance  parameters  p,  a  and  k  (in  some  applications,  these  might  be  known  a  priori). 
Once  the  reflectance  map  is  known,  we  must  then  solve  for  the  surface  orientation  along  the 
zero-crossing  contours,  assuming  that  p,  a  and  k  arc  constant  along  these  contours.  More  precisely, 
we  solve  the  following  problem.  Let  the  recorded  brightnesses  in  the  left  and  right  eyes  be 

Cl  =  #»[(1  —  «)m  +  au*] 

Cr  =  p[(l  —  a)m+  av*j 

respectively,  where 

m  =*  n  •  s 
u  (s  +  Vl)  •  W 

v/2  v/1  +  s  •  yL 
„  =  (s  +  vr)  ■  " 

s/2  VI  +  s  -  vR' 

We  will  first  solve  for  p,  a  and  k  over  some  portion  of  the  image,  and  then  use  these  values  to 
solve  for  the  surface  orientation  along  zero-crossing  contours  in  this  region  of  the  image. 

3.  Solving  for  the  Reflectance  Parameters 


The  task  of  determining  the  reflectance  parameters  would  be  considerably  eased  if  the  surface 
orientation  (p,  q)  were  already  known.  We  thus  restrict  our  attention  to  situations  in  which  this  is 
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true.  Consider  a  point  of  high  curs  attire  along  a  zero-crossing  contour,  vshich  in  the  limiting  case 
might  be  a  corner  where  two  zero-crossing  contours  intersect.  The  orientation  of  the  contour  in 
the  image  plane,  along  either  direction  of  the  contour,  is  different  (within  the  resolution  of  the 
image  grid).  A  directional  derivative  taken  along  one  direction  of  the  contour  yields  a  constraint 
on  the  relationship  between  p  and  q  at  this  point.  Since  there  arc  two  different  directions  along 
the  contour  at  such  a  point,  we  get  two  different  constraints  on  the  surface  orientation,  each  one 
corresponding  to  a  straight  line  in  gradient  space,  defining  a  linear  relationship  between  the  values 
of  p  and  q  at  the  point.  The  intersection  of  these  two  lines  uniquely  defines  the  surface  orientation 
at  the  comer. 

We  can  also  perform  this  computation  in  a  somewhat  non-standard  representation,  by  using 
unit  normals.  The  representation  we  choose  to  utilize  is  obtained  by  orthonormally  projecting 
the  Gaussian  hemi-sphere  of  all  visible  unit  normals  onto  the  plane  spanned  by  the  first  two 
components  of  the  unit  normal.  (Note  that  this  is  a  one-to-one,  onto  projection,  since  we  are 
only  interested  in  the  half  of  the  Gaussian  sphere  composed  of  the  unit  normals  corresponding 
to  visible  surfaces).  Thus,  rather  than  performing  computations  in  gradient  space,  we  instead  use 
the  components  of  the  unit  normals  lying  in  the  unit  disk.  In  this  ease,  each  directional  derivative 
constrains  the  first  two  components  of  the  unit  normal  to  lie  on  a  half  ellipse  in  the  unit  disk.  The 
intersection  of  two  such  half  ellipses  uniquely  defines  the  components  of  the  unit  surface  normal 
at  the  point. 

3.1.  Using  Irradiance  Derivatives 

We  now  consider  methods  for  determining  the  surface  reflectance  parameters  at  such  a  point. 
Initially  we  have  three  unknowns  and  only  two  equations.  To  uniquely  solve  for  the  unknowns, 
we  will  need  additional  information.  One  possibility  is  to  consider  using  derivatives  of  the  image 
irradiances  to  solve  for  the  parameters.  Taking  first  directional  derivatives  of  the  irradiance 
equations  introduces  three  new  variables,  p,,p„(=  <?,)  and  qy.  (Note  that  q,  is  assumed  to  be 
equal  to  py,  since  we  shall  assume  that  the  surfaces  are  at  least  twice  continuously  differentiable.) 
This  gives  us  a  total  of  8  unknown  variables,  the  three  listed,  as  well  as  p,  g,  p,  a  and  k.  At  a 
comer,  we  can  compute  the  values  of  p  and  q.  Furthermore,  second  directional  derivatives  along 
the  depth  contour  yield  two  additional  linear  constraints  on  pr,py  and  qv,  so  there  arc  only  4 
independent  variables.  We  also  have  four  equations,  namely: 

—■  —  P(1  —  a){mppx  4  m,pv}  4-  pak uk~l  {uppx  -f  u,p„} 

QP 

-T -  =  pO  -  <*MmpPv  +  "Mv)  +  p*kvk~'  {fcpPp  +  tM„} 

ay 

ap 

~  =  p(l  -  a){mppx  +  mqpy}  +  pakvk-'{vppz  +  vqpv) 

ap 

=  PO  -  “H"»pPv  +  mqqv)  4  pakvk~l{vppv  4  vqq„) 

ay 

where  the  subscripts  L  and  R  denote  the  left  and  right  images,  and  the  other  subscripts  denote 
partial  differentiation. 
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3.1  Solving  The  Equations 


We  now  seek  a  closed  form  so',  ion  to  this  set  of  equations;  that  is.  we  seek  expressions  for 
determining  p ,  a  and  ik  as  functions  of  the  image  irradiancc  derivatives,  the  illuminant  direction 
and  the  vergcnce  angle.  While  it  is  possible  to  obtain  such  a  closed  form  solution  (for  completeness, 
the  solution  is  given  in  the  appendix),  it  is  not  very  useful  for  computational  purposes,  since  it  is 
numerically  very  sensitive,  as  discussed  below. 

3.3.  Numerical  Stability 

While  the  appendix  describes  a  solution  to  the  problem  of  determining  the  parameters  of  the 
surface  reflectivity,  numerical  simulations  indicate  that  the  solution  is  not  numerically  stable.  ITiis 
is  illustrated  by  the  following  example,  used  to  determine  whether  the  image  values  &(*'/.  -  c») 
and  -§^{EL  —  Eh)  can  be  reliably  measured  from  an  image. 

Since  a  visual  system  is  not  an  infinite  resolution  device,  it  is  important  to  take  this  into 
account  when  simulating  the  solution  to  the  above  equations.  In  particular,  one  cannot  use 
analytic  expressions  for  terms  such  as  &(£i.  -  Eh)  but  rather,  one  must  use  finite  difference 
approximations  to  such  differential  expressions,  in  order  to  reflect  the  discrete  computational 
structure  of  the  vision  system.  Moreover,  the  range  of  values  for  £,  and  Eh  should  also  be  discrete, 
rather  than  continuous,  since  any  sensor  will  have  a  finite  resolution  of  brightness  values.  These 
considerations  turn  out  to  be  critical  in  evaluating  the  numerical  stability  of  the  computation. 
In  particular,  the  closed  form  solution  for  the  reflectance  parameters  given  in  die  appendix  in 
principle  provides  exact  values  for  those  parameters,  if  infinite  precision  data  is  available.  When 
discrete  approximations  arc  used,  however,  this  exact  solution  is  no  longer  possible,  as  is  illustrated 
by  the  following  example. 

A  synthetically  shaded  sphere  of  radius  160  pixels  with  a  total  irradiancc  range  of  8  bits  (or 
256  grey  levels)  was  constructed  using  equation  (1)  with  the  parameters  p  =  256,  a  =  0.8  and  with 
k  varying  from  1  to  10.  Note  that  this  implies  that  recorded  brightness  was  descretircd  to  a  finite 
precision,  thereby  restricting  the  resolution  of  measurements  of  spatial  change  in  the  difference 
between  the  image  brightnesses. 

The  ratio  of  interocular  separation  to  fixation  distance  was  set  at  6: 100.  The  total  number 
of  points  visible  to  both  eyes  was  counted,  and  the  percentage  of  such  points  for  which  the  finite 
differences  corresponding  to  the  above  partial  derivatives  exceeded  a  threshold  of  1  grey  level,  or 
1  part  in  256,  was  measured.  This  percentage  was  observed  to  be  small,  no  more  than  4%.  even 
for  large  values  of  k.  This  held  true  over  a  range  of  illuminant  directions.  Some  examples  are 
shown  in  Table  1.  (Similar  results  hold  for  other  values  of  a.) 


Table  I:  Sensitivity  in  computing  derivatives  of  image  differences 

The  table  lists  the  percentage  of  points  on  the  sphere,  visible  in  both  eyes,  for  which  the  finite 
difference  approximation  of  the  partial  derivative  of  the  change  in  image  intensity  exceeded  a  threshold  of 
one  part  in  256.  Here,  k  denotes  the  exponent  of  the  specular  component,  and  vis  is  the  percentage  of  all 
possible  points  on  the  sphere  which  are  visible  to  both  eves.  The  first  two  components  of  the  unit  source 
normal,  p.  and  q>.  are  given  by  p,  =  0,  ?,  —  0  in  the  left  columns,  and  ps  0.557,  =  0.239  in  the 

right,  and  the  vergence  angle  is  p  -  tan-1  j^. 

This  observation  suggests  that  the  system  of  equations  listed  in  die  appendix,  while  providing 
a  solution  in  the  ideal  case,  will  not  provide  a  numerically  stable  solution.  One  explanation  for 
this  can  be  seen  by  noting  that  any  partial  derivatives  of  the  irradiance  equations  will  involve 
partial  derivatives  of  die  components  of  the  surface  normal.  Unless  the  surface  has  a  significant 
curvature  at  diese  points,  the  partial  derivatives  of  the  image  irradiances  will  he  small  and  hence 
the  computation  of  the  surface  parameters  will  be  numerically  unstable.  For  example,  if  a  smaller 
sphere  had  been  used  in  the  above  example,  the  numerical  stability  of  the  solution  to  the  reflectance 
parameters  would  be  improved,  since  die  effective  curvature  of  the  sphere  would  be  increased. 

While  using  finite  differences  in  the  change  in  image  brightness  between  die  two  views  did 
nut  provide  detectable  values,  it  was  observed  that  the  simple  image  differences  frequently  did. 
Kxamplcs  are  shown  in  Table  2.  Here,  we  observe  that  even  for  strongly  specular  surfaces  (it  =  10) 
the  percentage  of  orientations  with  a  detectable  (>  1  grey-level)  difference  in  image  brightnesses 
is  over  55%.  Allowing  for  image  sensor  noise  (for  example,  image  differences  >  3)  still  leaves  an 
appreciable  percentage  of  orientations  with  measurable  differences. 


Table  2:  Sensitivity  in  computing  image  differences. 

The  table  lists  the  percentage  of  points  on  the  sphere,  visible  in  both  eyes  with  a  recorded  intensity 
of  at  least  5  parts  in  2S6.  for  which  (he  difference  in  image  intensity  exceeded  a  given  threshold.  Here,  k 
denotes  the  exponent  of  the  specular  component,  p  =  tan-1  and  p,  =  0.557  and  -  0.239.  The  first 
column  lists  the  percentage  of  points  on  the  sphere  which  are  visible  in  both  eyes  with  sufficient  brightness. 
The  remaining  three  columns  list  the  percentage  of  such  points  which  have  an  irradiance  difference  in  the 
two  eyes  of  the  amount  listed. 

This  suggests  that  while  we  may  not  be  able  to  rcliabiy  measure  the  spatial  derivatives  of 
the  difference  in  image  brightness,  due  to  the  discrete  nature  of  the  sensing  of  brightness  values 
and  the  discrete  approximations  to  the  spatial  derivatives,  we  may  be  able  to  rely  on  the  simple 
differences  in  image  brightness.  We  thus  turn  to  the  problem  of  finding  a  solution  for  die  surface 
reflectivity  parameters  from  the  image  brightnesses  alone. 

4.  A  Second  Solution 


We  first  observe  that  even  recording  the  irradianccs  at  a  corner  (defined  here  are  a  point  of 
high  curvature),  where  the  surface  orientation  is  known,  is  not  sufficient  to  solve  the  problem,  since 
we  have  two  equations  in  three  unknowns.  We  thus  require  an  additional  source  of  information 
and  will  assume  that  we  know  the  surface  orientation  at  two  distinct  points,  (we  will  show  later 
that  these  could  be  two  nearby  comers  or  simply  one  comer  and  a  nearby  point). 


At  each  comer,  we  note  that  the  following  expressions  hold  (where  i  denotes  the  image  point 
in  question) 


a. 


_ 2m,(fz.  —  £r)x _ 

(£ L  +  £«),(**?  -  vf)  -  (£c  -  £«),(«,*  +  ~  2m.) 

1 


1  + 


if*  ).«?-< eL ).«* ' 


-  +  -  «f)  -  (ft-  -  +  vf  -  2m,) 

fJ‘  2m.  (uf  -  vf) 

_  (gt  —  i  1 
u*  —  vk  a, 

This  gives  us  a  solution  for  p  and  a,  presuming  we  know  the  value  of  ft.  To  determine  fc.  we 
assume  that  the  surface  reflectivity  parameters  are  constant  over  the  region  of  the  image  spanned 
by  the  two  known  points.  In  this  case,  a,  =  a2  and  pi  =  p2  and  this  leads  to  the  following 
implicit  equations  for  k: 


(gr.  -  £ h),  _  rn3(ukCL  - 
{£l  —  £r)2  ~  rn,( -  vk£R)3 

and 

(£l  -  £r),  («*-«*), 

(fr.  -  £h)2  K-v*)2 

Hither  equation  can  be  solved  numerically  to  determine  the  value  of  k  and  hence  of  the  other 
reflectivity  parameters  a  and  p. 
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4.1.  Numerical  Stability 


While  the  equations  derived  in  section  3  were  observed  to  be  numerically  unstable,  the 
above  solution  is  much  more  sound.  This  can  be  observed  by  the  following  simulation.  Consider 
a  sphere  illuminated  by  a  point  source,  and  observed  by  two  sensors  whose  angle  of  vcrgence  is 
0  =  tan ~ 1  ( xlfo )  (roughly  equivalent  to  a  human  observer  fixating  on  an  object  1  meter  removed). 
1116  surface  of  the  sphere  can  be  given  any  surface  reflectance  properties  by  specifying  the 
parameters  k,p  and  a.  For  some  set  of  parameters,  we  compute  the  observed  brightness  value  in 
each  sensor  to  one  part  in  256.  We  then  apply  the  equations  derived  above  (assuming  initially 
that  the  value  of  k  is  known)  to  compute  the  parameters  a  and  p.  The  accuracy  of  the  computed 
values  is  compared  to  the  original  value  under  a  number  of  situations. 

Since  a  certain  amount  of  error  will  follow  from  the  truncation  of  brightness  accuracy,  we 
measured  the  error  in  computing  the  parameters  as  a  function  of  the  difference  in  observed 
brightness  values.  It  will  be  noted  that  we  can  trade  off  accuracy  of  computed  parameters  with  the 
percentage  of  orientations  for  which  such  a  computation  can  be  made  (see  t  able  3).  In  particular, 
as  we  increase  the  required  lower  bound  on  the  observed  differences  in  brightness,  we  decrease  the 
total  percentage  of  orientations  for  which  a  computation  can  be  made,  but  increase  the  percentage 
of  such  points  which  give  rise  to  computed  surface  parameters  within  a  particular  error  range, 
l  ablc  3  illustrates  some  examples. 
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Table  3:  Errors  in  computing  a  and  p. 
The  following  parameters  are  defined: 
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vis  is  the  percentage  of  viewed  orientations  on  a  sphere  for  which  the  difference  in  observed  brightnesses 
exceeds  some  threshold; 

a-t  is  the  percentage  of  vis  which  lie  within  a  range  ±1%  of  the  correct  value  of  a; 
p-i  is  the  percentage  of  vis  which  lie  within  a  range  ±e%  of  the  correct  value  of  p; 

(  is  the  percentage  error  bar  allowed  about  the  correct  value. 

In  this  case,  the  illuminanl  direction  is  given  by  p.  =  0.557,  =  0.239.  In  the  first  table,  the 

threshold  of  image  differences  was  set  at  3.  in  the  second  table,  the  threshold  was  5. 

To  compute  the  value  of  He,  we  must  solve  an  implicit  equation  in  k.  This  can  be  done 
numerically,  for  example,  by  Newton-Raphson.  While  computer  simulations  indicate  that  this 
compulation  is  very  robust  in  the  case  of  perfect  image  brightnesses,  the  computation  is  much 
less  stable  when  considering  truncated  brightness  values  (i.e.  if  infinite  precision  brightness  values 
are  used,  the  computation  returns  very  accurate  values  for  k ,  but  when  the  finite  resolution  of 
the  brightness  sensor  is  taken  into  account  the  accuracy  of  the  computation  drops  off  quickly). 
Table  4  indicates  an  example  of  the  error  associated  with  computing  k  numerically.  Ibis  table  was 
computed  by  choosing  two  surface  orientations  at  random,  such  that  the  difference  in  brightness 
exceed  a  threshold  of  3  grey  levels,  and  solving  the  implicit  equation  for  k  by  a  Newton-Raphson 
method. 

One  possible  method  for  improving  the  accuracy  of  computing  not  only  k,  but  also  p  and 
q,  is  to  consider  the  reflective  properties  of  the  surface  to  be  constant  over  some  extended  region 
of  the  image.  In  this  case,  a  number  of  points  may  be  sampled,  and  the  computed  values  for  the 
parameters  averaged.  Table  4  illustrates  the  reduction  in  error  associated  with  computing  k. 
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Table  4;  Errors  in  computing  k. 

Two  different  orientations  were  chosen  at  random,  such  that  the  difference  in  brightness  values  at 
each  point  exceeded  a  threshold  of  3  grey  levels.  The  percentage  of  such  pairs  of  points  with  percentage 
error  in  computing  k  within  the  given  range  are  indicated.  In  the  case  illustrated  below,  the  illuminanl 
direction  was  given  by  p.  =  0.557,  q,  =  0.239.  with  p  =  256,  a  =  0.5  and  k  =  5.  The  percentages  are 
shown  for  the  case  of  sampling  1.  3.  S,  and  10  pairs,  and  averaging  the  results. 


5.  Solving  for  the  Surface  Orientation 

Given  that  the  equations  derived  in  Section  4  can  be  used  to  determine  the  surface  reflectivity 
parameters  *,  p  and  a,  we  turn  to  the  final  stage,  which  is  to  compute  the  surface  orientation 
along  the  matched  zero-crossing  contours  of  the  stereo  algorithm. 


The  algorithm  described  below  is  essentially  a  modification  of  the  photometric  stereo  algorithm 
[Woodham,  1978,  1980,  1981:  Ikeuchi  and  Horn,  1979;  Silver,  1980J.  It  can  be  sketched  as  follows. 

(1)  Compute  the  reflectance  map  [Horn,  1977,  1981;  Horn  and  Sjobcrg,  1979]  associated 
with  the  computed  values  of  p,  a  and  k.  Since  R  is  a  function  of  the  unit  surface  normal, 
it  can  be  calculated  over  the  unit  disk,  with  coordinate  axis  given  by  the  first  and  second 
components  of  the  unit  normal. 

(2)  Use  the  measured  value  of  <f/  at  a  point  to  determine  a  contour  of  isobrightness  in  the 
unit  disk.  This  defines  the  set  of  possible  surface  orientations,  which  are  consistent  with  the 
value  of  £l- 

(3)  Use  the  measured  value  of  £h  at  the  corresponding  point  to  determine  a  second  contour 
of  isobrightness  in  the  unit  disk.  Ihc  intersections  of  the  two  contours  define  possible  surface 
orientations.  Note  that  in  genera!  there  will  be  at  least  two  such  points  of  intersection. 

(4)  Wc  can  disambiguate  the  possible  orientations  by  applying  a  third  constraint.  Take  a 
directional  derivative  of  depth  along  the  contour  through  the  point.  This  will  yield  a  linear 
constraint  between  p  and  q.  When  translated  to  the  unit  disk  (and  a  constraint  between  the 
components  of  the  unit  surface  normal),  this  constraint  becomes  a  half  ellipse  in  the  unit 
disk.  The  intersection  of  this  curve  with  the  two  isobrightness  contours  defines  the  correct 
surface  orientation  for  the  point  in  question. 

The  key  question  still  to  be  considered  'S  whether  this  algorithm  results  in  a  unique  solution. 
Wc  first  observe  that  the  isobrightness  contours  for  either  a  matte  or  a  specular  surface  consist  of 
a  scries  of  nested  ellipses  or  half-ellipses.  In  particular,  for  a  given  value  i?(n)  —  c,  the  ellipse  has 
center  in  the  unit  disk  of  (cp,,cqs)  where  ps  and  qs  are  the  first  and  second  components  of  the 
unit  source  vector  s.  The  minor  axis  of  the  ellipse  is  oriented  along  the  line  from  the  origin  to 
the  projection  of  s  and  the  major  axis  is  perpendicular  to  it. 

While  it  is  possible  to  derive  this  analytically,  it  more  easily  observed  by  noting  that 
for  a  Lambertian  coated  sphere,  with  source  located  along  the  same  direction  as  the  viewer, 
the  isobrightness  contours  arc  nested  circles.  For  an  off-axis  light  source,  the  corresponding 
isobrightness  contours  can  be  obtained  by  rotating  the  sphere  relative  to  the  viewer,  and  then 
projecting  the  previous  isobrightness  contours  onto  the  unit  disk,  clearly  resulting  in  a  series 
of  nested  ellipses.  A  similar  argument  holds  for  the  isobrightness  contours  of  a  purely  specular 
surface. 

When  considering  a  convex  combination  of  two  such  reflectance  maps,  the  isobrightness 
contours  become  somewhat  more  complex.  In  general,  however,  the  contours  are  slightly  distorted 
ellipses  (sec  Figure  3).  As  a  consequence,  the  intersection  of  isobrightness  contours  from  the  two 
eyes  will  in  general  result  in  a  pair  of  points  in  the  unit  disk.  (This  excludes  the  degenerate  case 
of  a  pure  matte  surface,  in  which  case  the  intersection  will  be  a  complete  ellipse.)  We  must  now 
determine  the  circumstances  under  which  the  third  constraint  will  disambiguate  this  pair  of  points. 

Suppose  that  the  zero-crossing  contour  at  the  point  of  interest  has  a  tangent  in  the  image 
plane  whose  angle  with  respect  to  the  z  axis  is  7,  and  whose  derivative  of  depth  along  the  contour 


Figare  3.  Examples  of  isobrightness  contours  on  the  unit  disk,  for  different  reflectance  maps.  In  all 
four  cases,  p,  =  0.7,  q.  —  0.3.  For  cases  (a),  (b),  and  (c).  the  exponent  k  =  to,  and  the  parameter  a  is 
0.1, 0.5  and  0.9  respectively.  In  case  (d),  the  exponent  is  k  =  5  and  a  =  0.5. 


Figure  4.  Constraints  on  the  unit  surface  normal,  as  a  function  of  the  component  of  the  surface 
gradient  along  a  contour.  Kach  contour  delineates  the  possible  values  of  the  unit  surface  normal,  for  a 
different  value  of  the  surface  gradient  along  the  zero-crossing  contour.  In  each  of  the  contours  illustrated 
here,  the  parameter  7  was  negative.  The  straight  line  corresponds  to  the  case  in  which  the  depth  derivative 
along  the  zero-crossing  contour  w<js  zero. _ _ 

is  given  by  the  value  d.  This  derivative  determines  the  component  of  the  surface  normal  along  the 
contour.  When  translated  into  a  constraint  on  the  components  of  the  unit  surface  normal,  the  set 
of  possible  unit  normals  which  arc  consistent  with  such  a  component  along  the  contour  are  given 
by  the  parametric  equations  (assuming  7  is  non-zero) 


n,(i)  = 


na(t)  = 


t  sin  7 

(»in2  7  d2  —  2dt  co»7  -f  t2)* 
d  —  t  cos  7 

(*in27  -f-  d2  —  2dt  cos 7  -f- 12)^ 


These  parametric  equations  define  a  half  ellipse  in  the  unit  disk  with  center  at  the  origin  and  with 
major  axis  the  line  through  the  origin  with  slope  —  cot  7  of  length  2.  If  7  is  positive,  the  half 
ellipse  is  that  which  connects  the  two  endpoints  of  the  major  axis  and  whose  half  of  the  minor 
axis  lies  the  the  positive  n2  half  disk,  if  7  is  negative,  the  half  ellipse  lies  with  its  half  of  the 
minor  axis  in  the  negative  n2  half  disk  (see  Figure  4). 


We  wish  to  determine  (he  conditions  in  which  there  is  not  a  unique  solution,  which 
corresponds  to  determining  the  conditions  under  which  the  half  ellipse  constraint  passes  through 
both  points  arising  from  the  intersection  of  the  isobrightness  contours.  Note  that  the  family  of  all 
possible  depth  constraints  can  be  represented  by  the  pairs  (7,6)  where  7  ranges  from  —r  to  * 
and  b  is  the  length  of  the  minor  axis  of  the  half  ellipse.  Consider  a  point  (po,fo)  lying  in  the  unit 
disk.  L^t 

(  —  sgn(f0)Un-1  — . 

Po 
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Then  the  valid  range  of  possible  7  is  from  £  —  J  to  £  +  $,  modulo  2k.  Thus,  the  set  of  possible 
depth  constraints  which  pass  through  this  point  can  be  specified  by  the  graph  given  in  Figure  5. 
which  has  a  range  of  *  in  -7  and  a  maximum  in  6  of  \'p\  +  ?jj.  For  a  second,  distinct  point  in 
the  unit  disk,  a  similar  mapping  is  obtained,  which  intersects  the  first  at  exactly  one  point.  This 
implies  that  there  is  exactly  one  ellipse,  with  orientation  and  eccentricity  given  by  this  intersection, 
which  passes  through  both  points  in  the  unit  disk.  Ihis  in  turn  implies  that  the  depth  constraint 
will  almost  always  uniquely  specify  the  surface  orientation  of  the  zero-crossing  point,  since  only 
a  particular  alignment  of  the  light  source  and  the  depth  contour  will  result  in  an  ambiguity,  and 
in  general,  the  physical  processes  which  give  rise  to  the  depth  contour  constraint  arc  independent 
of  the  light  source.  In  the  case  of  continuous  value  ranges  for  7  and  d.  the  set  of  conditions 
in  which  a  unique  answer  is  not  possible  will,  in  fact,  have  measure  zero.  In  the  more  practical 
case  of  a  discrete  range  of  values  for  7  and  d,  assume  that  the  orientation  7  can  take  on  one 
of  n  different  values,  while  the  component  of  the  depth  gradient  along  die  contour  d  can  take 
on  one  of  m  different  values.  Then  assuming  that  the  possible  values  for  7  and  d  arc  uniformly 
distributed,  and  assuming  that  the  constraints  due  to  the  reflectance  map  isobrightness  contours 
arc  independent  of  the  constraint  due  to  the  depth  derivative  implies  that  the  probability  of  an 
ambiguous  computation  of  surface  orientation  is  Clearly,  if  we  have  a  reasonable  resolution 
for  7  and  d ,  this  probability  will  be  very  small. 

6.  Discussion 

We  note  that  it  is  possible  to  solve  for  the  reflectivity  parameters  without  requiring  two 
comers.  Suppose  instead  that  we  have  one  comer,  and  that  the  surface  is  at  least  twice  continuously 
differentiable.  Then  by  noting  that  integration  of  surface  orientation  around  a  closed  loop  must 
yield  zero  (note  that  a  similar  constraint  has  been  used  by  Ikcuchi  and  Horn  [1981])  we  can 
determine  the  surface  orientation  at  a  second  point,  and  the  analysis  of  Section  4  still  holds. 

We  also  note  that  the  numerical  properties  of  the  algorithm  sketched  here  probably  preclude 
its  use  by  the  human  visual  system.  In  particular,  to  have  reasonable  numerical  stability  over  a 
large  range  of  surface  orientations,  a  large  vergcnce  angle  is  required  (on  the  order  of  3°).  For  the 
human  system,  this  is  roughly  equivalent  to  viewing  objects  from  a  distance  of  1  meter  or  less. 
In  cases  of  machine  vision,  where  control  of  the  positioning  of  the  sensors  is  possible  (i.e.  we  arc 
not  restricted  to  interocular  separations  of  6  cm)  large  vergence  angles  are  more  feasible. 

7.  Error  Analysis 

As  well  as  performing  simulations  of  the  computation  of  surface  parameters,  we  can  also 
make  an  estimate  of  the  error  associated  with  computing  the  parameters.  Ihis  analysis  will  also 
indicate  a  means  for  determining  those  situations  in  which  the  error  is  small. 

We  again  consider  the  situation  in  which  the  parameter  k  is  known  exactly,  as  well  as  the 
values  of  p  and  q  and  hence  of  m,  u  and  v.  We  let 


r 


* 


Figure  5.  Determining  the  uniqueness  of  the  solution.  Parts  (a)  and  (b)  represent  constraints  between 
the  parameters  6  and  7,  for  different  values  of  {.  Wien  a  second  point  is  added  (c),  the  intersection  of 
the  two  contours  defines  the  only  possible  half  ellipse  in  the  unit  disk  which  passes  through  both  points 
in  the  unit  disk  and  has  appropriate  major  and  minor  axes.  Thus,  given  two  isobrightness  contours  in  the 
unit  disk,  from  the  two  views  of  a  point  on  a  zero-crossing  contour,  there  is  exactly  one  depth  contraim 
0  which  will  not  disambiguate  these  points. 
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a,„  =  measured  value  for  a 
pm  =  measured  value  for  p 
a0  —  actual  value  for  a 
Po  =  actual  value  for  p 
it  =  error  in  £L 
Sr  =  error  in  £R. 

The  goal  is  to  determine  estimates  for  the  ratio  of  the  measured  values  to  the  actual  values,  that 
is,  to  derive  expressions  for 


- 

£»0 


By  expansion,  we  obtain  the  following  expression 


Thus, 


Similarly, 


Pm  = 


(e 


El  —  Eh  +  &l  —  Sr' 


: ik  —  «* 


x 

(1  4-  Sl  —  Sr  \  Po 
El  —  Er)  t. 


_ MEl  —  Er  +  6l  —  Sr) _ 

°m  (t*fc  —  "»Xf«  +  Sr)  —  (vk  —  m)(£L  +  SL) 

f. *  ,  -  6r  _ («*  ~  ~  K  ~  m)6L 

°\  El-Er)[  («*  -  mXfR  +  SR)  -  (**  -  m)(£L  +  «,.) 

so  that 

f.  ,  6l-6r\( _ («*  ~  ~  (»*  ~  ™)Sl  \ 

V  El  —  Er)\  {uk  —  6r)  —  (vfc  —  m)(£L  +  6L)J 


The  goal  is  to  determine  what  circumstances  will  result  in  small  errors  in  the  computed 
parameters  (i.e.  under  what  conditions  will  ep  ss  1  and  ta  «  1). 


Algebraic  substitution  yields  the  following  alternative  forms: 

,  ,  Sl  —  Sr  1  Sr£l  —  6lEr 
P  El  —  Er  p  —  Er) 

.  Sl-Sr _ Sl  Sr  -f-  ^(ElSr  —  ErSl)  A 

\  El  Er  /  y  El  —  Er  6l  —  Sr  4-  ElSr  —  Er6l)  J 

We  can  see  that  the  only  circumstances  in  which  the  error  will  be  large  arc  when  \El  —  Er\  is 
not  significantly  greater  than  | SL  —  6«|.  Since  in  general,  | Sl  —  $k|  will  be  no  more  than  1.  by 
restricting  our  attention  to  portions  of  the  image  in  which  )EL  —  En\  is  larger  than  some  threshold, 
wc  can  reduce  the  probable  error  in  our  estimation  of  the  surface  reflectivity  parameters. 
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8.  Summary 


Wc  have  demonstrated  that  shading  information  can  be  combined  with  stereo  data  along 
feature  point  contours  in  order  to  determine  the  surface  orientation  along  those  contours,  under 
the  following  conditions. 

(1)  The  surface  is  illuminated  by  a  distant  point  source,  whose  direction  relative  to  the  viewer 
is  known. 

(2)  'Hie  surface  is  viewed  from  two  different  positions,  and  the  angle  of  vcrgcncc  of  the  two 
views  is  known.  For  numerical  stability,  this  angle  should  be  on  the  order  of  3  degrees  or 
larger. 

(3)  The  surface  is  at  least  twice  continuously  differentiable. 

(4)  Ihe  reflectance  map  of  the  surface  can  be  represented  in  the  form 

R(n)  —  (1  —  a)(n  s)  4-  a{n  h}* 

where  s  is  a  unit  vector  in  the  direction  of  the  light  source,  h  is  a  unit  vector  in  the  direction 
of  the  off  specular  angle. 


and  n  is  a  unit  surface  normal.  The  variables  a,p  and  k  are  unknown,  and  are  to  be 
computed  from  the  irradiancc  data.  It  is  assumed  that  a  is  non-zero,  so  that  the  surface  is 
not  perfectly  matte. 

Under  these  conditions,  we  have  shown  that  using  the  brightness  values  recorded  at 
corresponding  points  in  the  two  images  (determined  by  the  stereo  correspondences)  we  can 
compute  the  reflectivity  parameters  for  the  surface,  and  using  in  addition  the  component  of  the 
surface  gradient  along  the  stereo  contour,  we  can  apply  an  algorithm  similar  to  photometric  stereo 
to  compute  the  surface  orientation  along  the  stereo  contours. 

The  numerical  stability  of  the  computation  was  also  explored  and  methods  for  reducing  the 
possible  error  in  the  computation  were  suggested.  The  key  observation  concerning  the  computation 
was  that  while  accurate  solutions  were  possible  given  high  resolution  input  (e.g.  the  image 
brightnesses  were  known  to  extremely  high  accuracy),  when  the  discrete  nature  of  the  brightness 
sensistivity  and  spatial  sensitivity  of  the  imaging  devices  was  taken  into  account,  the  numerical 
stability  of  the  computations  degraded  rapidly.  This  suggest  that  in  its  current  form,  computations 
such  as  those  outlined  here  may  be  of  limited  practical  use,  until  methods  for  improving  the 
numerical  stability  are  derived. 
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10.  Appendix 


The  solution  for  the  reflectivity  parameters,  using  derivatives  of  the  image  irradianccs  is 
summarized  as  follows.  The  solution  to  k  is  given  by 

In  {(o2  4  mpa3)u,  —  (ai  4  m,a3)t<p}  —  In  {(a2  —  mpa3)vg  —  (ai  —  mqa3)v p} 

In  v  —  lnu 

the  albedo  p  is  given  by 


2  (  *  ^  ±  Ev)  -  K  +  uD\*  ~  iO  ~ 

p  (/ifnp  wm,)4  +  [(1  -  i/)m„  +  (1  -  Om,]* 

and  the  specular  scalar  a  is  given  by 


where 


a  — 


i  * 

p  e 


a,  ==  jwipfAjEj,  —  AyE,)  4  m,(/iAyEv  —  wAxE*)  4  A  yE*m,(ji£  —  uv) 
a2  =  umq( A*Ey  —  AyE*)  4  mp(/iA„Ey  —  wA*Ej)  4  &zY,ymp(pZ  —  u/i/) 
a3  =  uA2z  4  (uiu  —  ni) A*Ay  — 


where 


and  where 


•  =  (A,  +  v&M  ~  (Ay  4  t*z)B 

♦  =  nAAy  —  ljBAx 

8  =  h-A?  (wy  —  p£)AB  —  u>B* 


dEt 

dx 

d€L 

dy 

dCL 


dx 

dy 

dEn 


dx  dx 
dEL.dER 
dy  dy 


v 

A  =  *(«,«*-> -v,*-) 
B 
C 
D 


k(uquk  1  -  v,v*  ‘) 
k(upuk~'  4  WpV*-1)  —  2mp 

k[uq u*— 1  4  vX“‘)  —  2m, 


and  where  the  second  directional  derivatives  along  the  two  directions  of  the  zero-crossing  contour, 
whose  values  are  d\  and  d2  with  associated  directions  71  and  72,  determine  the  constants 
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d,  sin2  72  —  d-^ain2  7, 

sin  (7j  -f  72)  sin  (72  —  7,) 

ds  cos2  7!  —  d\  cos2  72 

sin  (7,  ~j-  72)  sin  (72  —  71) 

2  sin  7 1  sin  72 
v  =  -7—. - r 

8111(7,  4-  72) 

^  _  2  cos  7 1  cos  72 
sin  (7,  +72) ' 

Finally,  the  surface  orientation  at  this  point  is  given  by  the  three  parameters 

Jly  ~  (A,  +  »&V)A  -  (Ay  -F  i&x)B 

fzi  =  H  —  vfzy 

fyy  ~  u  Cfzy 

10.1.  Discussion  of  the  Equations 

Several  comments  on  the  form  and  derivation  of  the  equations  arc  appropriate.  First,  we  note 
that  in  deriving  the  above  expression  for  a,  the  form  obtained  is  valid  only  in  the  case  of  both 
A, +  A,  and  A*  —  A„  non-zero.  This  provides  an  interesting  side  effect,  since  these  expressions 
both  vanish  if  and  only  if  the  surface  material  has  only  a  matte  component,  and  no  (measurable) 
specular  component,  or  the  surface  is  a  plane.  If  the  surface  material  has  no  specular  component, 
o  =  0  and  this  is  exactly  the  expression  obtained  by  substitution  of  A,  =  0  and  Ay  =  0  into  the 
expression  obtained  for  a. 

Second,  we  note  that,  subject  to  possible  numerical  errors,  the  solution  for  the  surface 
reflectance  parameters  is  unique. 

Third,  we  note  that  it  is  important  to  consider  the  measurability  of  the  image  knowns.  In 
other  words,  can  we  extract  the  required  measurements  for  the  images  in  a  reliable  manner?  We 
first  note  that  the  zero-crossing  contours,  along  which  the  stereo  depth  information  is  known,  will 
frequently  outline  discontinuities  in  the  surface  or  in  the  surface  material,  for  example,  occluding 
edges  of  objects  or  changes  in  surface  albedo.  Since  we  cannot  measure  infinitesimal  derivatives 
in  the  image,  does  this  preclude  our  ability  to  measure  derivatives  of  the  images?  The  answer 
is  no.  since  we  can  again  rely  on  directional  derivatives.  For  example,  we  need  to  measure 
—  i a)  and  ^[£l  —  £/<)■  This  can  be  done  without  crossing  a  discontinuity  in  the  image 
irradianccs  by  considering  the  two  directional  derivatives  along  the  contour  obtained  at  a  comer  in 
the  zero-crossing  contour.  At  this  point,  we  have  two  different  directional  derivatives,  and  simple 
algebra  allows  us  to  then  determine  the  parametric  partial  derivatives. 
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